Routes to spatiotemporal chaos in the rheology of nematogenic fluids 
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With a view to understanding the "rheochaos" observed in recent experiments in a variety of 
orientable fluids, we study numerically the equations of motion of the spatiotemporal evolution of 
the traceless symmetric order parameter of a sheared nematogenic fluid. In particular we establish, 
by decisive numerical tests, that the irregular oscillatory behavior seen in a region of parameter 
space where the nematic is not stably flow-aligning is in fact spatiotemporal chaos. We outline 
the dynamical phase diagram of the model and study the route to the chaotic state. We flnd that 
spatiotemporal chaos in this system sets in via a regime of spatiotemporal intermittency, with a 
power-law distribution of the widths of laminar regions, as in H. Chate and P. Manneville, Phys. 
Rev. Lett. 58, 112 (1987). Further, the evolution of the histogram of band sizes shows a growing 
length-scale as one moves from the chaotic towards the flow aligned phase. Finally we suggest 
possible experiments which can observe the intriguing behaviors discussed here. 

PACS numbers: 61.30.-v,95.10.Fh,47.50.+d 



I. INTRODUCTION 



The intriguing rheological behavior of solutions of en- 
tangled wormlike micelles has been the subject of a large 
number of experimental and theoretical studies in re- 
cent years 0, These long, semiflexible cylindrical 
objects, whose length distribution is not fixed by chem- 
ical synthesis and can vary reversibly when subjected 
to changes in temperature, concentration, salinity and 
flow, have radii ^ 20-25 A, persistence lengths ~ 150 
A and average lengths upto several microns. Like poly- 
mers they entangle above a critical concentration and 
show pronounced viscoelastic effects. However, unlike 
covalently bonded polymers, these "living polymers" can 
break and recombine reversibly in solutions, with pro- 
found consequences for stress relaxation and rheology 
in the form of shear banding 0, 0, Q, and rheological 
chaos|i|li,|i[li[ll[l3. Measurements mill report 
monoexponential relaxation of the viscoelastic response 
in accordance with the Maxwell model of viscoelastic- 
ity. However, for wormlike micelles of CTAT0, ^3 
concentration 1.35 wt.%, the fit to the Maxwell model 
is very poor, and the Cole-Cole plot deviates from the 
semi-circular behavior expected in Maxwellian systems 
and shows an upturn at high frequencies. This deviation 
from Maxwellian behavior is possibly due to the compa- 
rable values of timescales associated with reptation {Trep) 
and reversible scission (rb) in this system unlike in other 
wormlike micellar systems where the differences in the 
time scales << Trep lead to a 'motional averaging' 
effect. Further, in the concentrated regime, when the 
mesh size of the entangled micellar network is shorter 
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than the persistence length of the micelles, orientational 
correlations begin to appear Q. In fact the nature of 
viscoelastic response and the development of long-range 
orientational order at high concentration play an impor- 
tant role in the non-linear rheology of wormlike micelles, 
in particular in shear banding transition and rheochaotic 
behavior 113. 

In this paper we explore th e dy namical phase diagram 
of the model studied in [lol [lTI |. with emphasis on the 
route to spatiotemporal chaos. Our primary finding is 
summarized in Fig. ^ which shows that this route is 
characterized by spatiotemporal intermittency. Before 
presenting our results in more detail, we cover some nec- 
essary background material. 

The application of large stresses and strains on worm- 
like micellar solutions can result in a variety of complex 
rheological behavior. Many dilute solutions of worm- 
like micelles exhibit a dramatic shear thickening behav- 
ior when sheared above a certain threshold rate, often 
followed by the onset of a fiow instabilitv [l6lll7l Hq . Ex- 
periments have observed shear-banded flow in wormlike 
micellar solutions with formation of bands or slip layers 
of different microstructures havin g ve ry different rheolog- 
ical properties [mill US miiilll. The shear banding 
transition is a transition between a homogeneous and an 
inhomogeneous state of fiow, the latter being character- 
ized by a separation of the fiuid into macroscopic domains 
or bands of high and low shear rates. It is associated 
with a stress plateau (above a certain critical shear rate 
7c where the shear stress a versus shear rate 7 curve is a 
plateau) in the nonlinear mechanical response. 

More recently, rheological chaos or "rheochaos" has 
been observed in experiments studying the nonlinear 
rheology of dilute entangled solutions of wormlike mi- 
celles formed by a surfactant C TAT H Hill 13. Under 
controlled shear rate conditions in the plateau regime, 
the shear stress and the first normal stress difference 
show oscillatory and more complicated, irregular time- 



2 




FIG. 1: (Color online) Space-time plots (spatial variation 
along abscissa and time along ordinate) of the shear stress 
for 7 = 4.0 and (a) = 1.11 (time-periodic, spatially homo- 
geneous), (b) to (d), Xk = 1.12,1.13,1.15 (spatiotemporally 
intermittent), (e) Xk = 1.22 (spatiotemporally chaotic), (f) 
and (g), Xk = 1.25 and 1.27 (chaotic to aligning) and (h) 
Xk — 1.28 (aligned) (colormap used: black(low shear stress) 
red yellow(high shear stress). Slices taken from system 
of size L = 5000. 



dependence. Analysis of the measured time series shows 
the existence of a positive Lyapunov exponent and a fi- 
nite non-integer correlation dimension characteristic of 
deterministic chaos. 

Occurrence of sustained oscillations often of an irregu- 
lar nature have also been reported in some other experi- 
ments on complex fluids in shear flow. Roux et a/.|25ll26l| 
have observed sustained oscillations of the viscosity near 
the non-equilibrium, layering transition to the "onion" 
state in a lyotropic lamellar system consisting of close 



compact assembly of soft elastic spheres psl l26l 127^. It 

has been conjectured that the presence of oscillations in 
the viscosity is due to structural changes in the fluid, 
arising out of a competition between an ordering mecha- 
nism that is driven by stress and a slow textural evolution 
which destroys the stress induced ordered state. "Elastic 
Turbulence" in highly elastic polymer solutions "28] and 
"Director Turbulence" in nematic liquid crystals in shear 
flow [2^ Is^l are two other examples of highly irregular 
low- Reynolds-number flows in complex fluids. Both these 
phenomena are characterized by temporal fluctuations 
and spatial disorder. Also worth noting is the observa- 
tion by Ramamohan et al. jsj] of rheochaos in numerical 
studies of sheared hard-sphere Stokesian suspensions. 

Many complex fluids have nonlinear rheological con- 
stitutive equations that cannot sustain a homogeneous 
steady flow. This material instability occurs when the 
stress vs. strain rate curve is non-monotonic in nature, 
admitting multiple strain rates 7 at a common stress a. 
Particularly for shear flow, it has been shown [s^l that ho- 
mogeneous flow is linearly unstable in a region where the 
incremental shear viscosity is negative, i.e., da/dj < 0. 
The system then undergoes a separation into two co- 
existing macroscopic shear bands at different shear rates 
arranged so as to match the total imposed shear gradient. 
Systems where the dynamic variables u or 7 are coupled 
to microstructural quantities may admit many other pos- 
sibilities - the flow may never be rendered steady in time, 
or it may become spatially inhomogeneous even erratic 
or both. Fielding and Olmsted study one such scenario!^ 
in the context of shear thinning wormlike micelles where 
the flow is coupled to the mean micellar length. 

Significantly, Grosso et al. 33] and Rienacker et al. 
[33. 135I I , find temporal rheochaos in the dynamics of the 
passively advected alignment tensor alone. They study 
the well-established equations of hydrodynamics for a ne- 
matic order parameter, with material constants corre- 
sponding to a situation where stable flow alignment is 
impossible. They consider only spatially homogeneous 
states 36], i.e., they study a set of ordinary differential 
equations for the independent components of the nematic 
order parameter, evolving in the presence of an imposed 
plane shear flow. They are thus not in a position to 
explore the implications of the observed chaos for shear- 
banding. 

Other theoretical approaches aimed at explaining the 
rheological chaotic oscillations in a wormlike micellar 
fluid include those by Gates et aL@. In the shear thick- 
ening regime Gates et al.^ propose a simple phenomeno- 
logical model for a fluid with memory and an under- 
lying tendency to form shear-banded flows, with only 
one degree of freedom - the shear stress. Recently Ara- 
dian and Gates 37, 33 have studied a spatially inhomo- 
geneous extension of this model, with spatial variation 
in the vorticity direction. Working at a constant aver- 
age stress ((t), they observe a rich spatiotemporal dy- 
namics, mainly seen in what they call "flip flop shear 
bands" - a low and a high unstable shear band sepa- 
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rated by an interface and periodically flipping into one 
another. For a certain choice of parameters they observe 
irregular time-varying behavior, including spatiotempo- 
ral rheochaos. As in our work, the key nonlinearities in 
|37l Is^ arise from nonlinearities in the constitutive re- 
lation, not from the incrtial nonlinearities familiar from 
Navier-Stokes turbulence. An important result l3^ . l3^ 
is that they are able to find complex flow behavior even 
when the stress vs. shear-rate curve is monotonic. In 
addition, they find rheochaos even in a few-mode trun- 
cation where well-defined shear bands cannot arise. 

In this paper we study a minimal model to explain the 
complex dynamics of orientable fluids, such as wormlike 
micelles subjected to shear flow. We show that the basic 
mechanism underlying such complex dynamical behavior 
can be understood by analyzing the relaxation equations 
of the alignment tensor of a nematogenic fluid, the under- 
lying idea being that wormlike micelles being elongated 
objects will have, especially when overlap is significant, 
a strong tendency to align in the presence of shear. We 
study equations of motion of the nematic order parame- 
ter in the passive advection approximation i.e., ignoring 
the effect of order parameter stresses on the flow pro- 
file which we take to be plane Couette, incorporating 
spatial variation of the order parameter. We calculate 
experimentally relevant quantities, e.g, the shear stress 
and the first normal stress difference, and show that, 
in a region of shear rates, the evolution of the stresses 
is spatiotcmporally chaotic. Further, in this region the 
fiuid is not homogeneously sheared but shows "dynamic 
shear banding" (banded flow with temporal evolution of 
shear bands). A careful analysis of the space-time plots 
of the shear stress show the presence of a large number of 
length scales in the chaotic region of the phase space of 
which only a few dominant ones are selected as one ap- 
proaches the boundary of the aligned phase. Finally we 
explore the routes to the spatiotcmporally chaotic state. 
The transition from a regular state (either temporally 
periodic and spatially homogeneous, or spatiotcmporally 
periodic) to a spatiotcmporally chaotic one occurs via a 
series of spatiotcmporally intermittent states. By calcu- 
lating the dynamic structure factor of the shear stress, 
and the distribution of the sizes of laminar domains, we 
can distinguish this intermittent regime from the spa- 
tiotcmporally chaotic and the regular states occurring in 
this model. Finally, we present a nonequilibrium phase 
diagram showing regions where spatiotcmporally regular, 
intermittent and chaotic phases are found. 

The paper is organized as follows. In the next section 
we introduce the model and describe in detail the spa- 
tiotemporal chaos that we observe, along with the routes 
to chaos. We then conclude with a summary and discus- 
sions of our results. Our main results on the spatioteni- 
poral nature of rheochaos have appeared in an earlier, 
shorter article(lQ|. 



II. SPATIOTEMPORAL RHEOLOGICAL 
OSCILLATIONS AND CHAOTIC DYNAMICS IN 
A NEMATOGENIC FLUID: 

A. Model and Methods 

Traditionally, complex rheological behaviors such as 
plateau in the stress vs shear-rate curveQ, shear- 
banding|3, 0, and "spurt" ^ have been understood 
through phenomenological models for the dynamics of 
the stress such as the Johnson-Segalman fJS);40l. |4]| 
model, which produce non-monotonic constitutive rela- 
tions. In such equations the stress evolves by relaxation 
or by coupling to the velocity gradient. For example in 
the JS model the non-Newtonian part of the shear stress 
(T evolves according to 

-^+u-'VcT+(T[n—aK] + [fl — aK]'^cr = 2fiK — Ta~^cr (1) 

with a stress relaxation time tq, an elastic modulus fj, 
and a parameter 'a' (called the slip parameter) con- 
trolling the non-affine deformation, u is the hydrody- 
namic velocity field and and k are the antisymmetric 
and symmetric parts of the rate-of-deformation tensor. 
A useful point of view, and one that unifies such phe- 
nomenological descriptions with dynamical models of or- 
dering phenomena in condensed matter physics, is that 
such equations of motion for the stress are not funda- 
mental but are derived from the underlying dynamics 
of an alignment tensor or local nematic order param- 
eter Q. Equations of motion for the latter are well- 
establishedgl EHillillillilliililli^in terms of 
microscopic mechanics (Poisson brackets) and local ther- 
modynamics, and naturally include both relaxation and 
fiow-coupling terms of essentially the sort seen, e.g., in 
the JS model. The contribution of the order parameter 
Q to the stress tensor is also unambiguous within such a 
framework, once the free-energy functional F[Q] govern- 
ing Q is specified. This approach is particularly appropri- 
ate when the system in question contains orientable en- 
tities, such as the elongated micelles of the experiments 
of 6, 7, 15, 24]. Thus, not worrying about properties 
specific to a wormlike micelle, e.g., the breakage and re- 
combination of individual micelles, one can attempt to 
understand the properties of the wormlike-micelle solu- 
tion by treating it as an orientable fluid and analyzing 
the equations of motion of the nematic order parameter. 
While properties specific to living polymers might play 
an important role in their rheological behavior, the gen- 
erality of our order parameter description encourages us 
to think that we have captured an essential ingredient for 
rheochaos. As we shall see, this approach leads naturally 
to terms nonlinear in the stress, absent in the usual JS 
equations of motion, which lead ultimately to the chaos 
with which this paper is concerned. Refs. [13, 113 found 
it necessary to modify the Johnson-Segalman equation 
by including terms nonlinear in the stress in order to 
produce chaos. 
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We now discuss the relaxation equation of the align- 
ment tensor characterizing the molecular orientation of 
a nematic liquid crystal in shear flow. These e quat ions 
were derived by various groups Qi llsl l47l IIM 

I49I ISQj . using different formalisms resulting in broadly 
similar though not in all cases identical equations of mo- 
tion. We work with the equations of'47', 'Sl'l, so as to 
make contact with the recent studies 52] of purely tem- 
poral chaos in the spatially homogeneous dynamics of 
nematic liquid crystals in flow. These authors have ex- 
tended their analysis 0, 113 to include biaxially ordered 
steady and transient states. Their work has revealed a 
transition from a kayaking-tumbling motion to a chaotic 
one via a sequence of tumbling and wagging states. Both 
intermittency and period doubling routes to chaos have 
been found. 

A nematogenic fluid is comprised of orientable objects, 
such as rods or discs, with the orientation of the ith par- 
ticle denoted by the unit vector Ui. In the nematic phase 
there is an average preferred direction of these molecules, 
which distinguishes it from the isotropic phase where 
there is no such preferred direction. The order parame- 
ter that measures such apolar anisotropy is the traceless 
symmetric "alignment tensor" or nematic order parame- 
ter 

1 ^ 1 
QMr) - - 5]((^.V^^ - -5^p))S{r - r^). (2) 

i=l 

built from the second moment of the orientational distri- 
bution function. By construction, it is invariant under 
Oi —v>i and vanishes when the Ui are isotropically 
distributed. 

Since nematic fluids possess long range directional or- 
der, presence of spatial inhomogcneitics would result in 
deformations of the director field and hence cost elastic 
energy. In general when the variable describing an or- 
dered phase is varying in space, the free energy density 
will have terms quadratic in VQ. Static mechanical equi- 
librium for Q corresponds to extremising the Landau-de- 
Gennes free-energy functional 

F[Q] ^ J d'xl^Q : Q - ^BiQ • Q) : Q + ^(Q : Q)^ 
+ yVQ; VQ + yV-Q- V-Q], (3) 

with phenomenological parameters A, B and C governing 
the bulk free-energy difference between isotropic and ne- 
matic phases, and Fi and r2 related to the Frank elastic 
constants of the nematic phase. A macroscopically ori- 
ented nematic with axis n has Q = s{nn — 1/3) (where I 
is the unit tensor) which defines the conventional (scalar) 
nematic order parameter s. For A small enough but pos- 
itive, F has minima at s = and at s = sq ^ 0. The 
minimum at s = Sq is lower than the one at s = when 
A < A^ = 2_B^/9C which corresponds to the (mean- field) 
isotropic-nematic transition. The functional F plays a 



key role in the dynamics of Q as well. The equation of 
motion for Q is 

— -+U-VQ = T~-^G + (aoK -f ai/« • Q)sT 

at 

+ rj Q Q rj, (4) 

where the subscript ST denotes symmetrization and 
trace-removal, u is the hydrodynamic velocity field, 
K = (l/2)[Vu-f-(Vu)'^] and = (l/2)[Vu- (Vu)^] the 
shear-rate and vorticity tensors respectively. The flow ge- 
ometry imposed is plane Couette with velocity u = ^yx 
in the x direction, gradient in the y direction and vortic- 
ity in the z direction, r is a bare relaxation time, ao and 
ai are parameters related to flow alignment, originating 
in molecular shapes. G, the molecular field conjugate to 
Q, is given by 

G = -(5F/5Q)sT = -[AQ-x^B(Q-Q)sT + CQQ:Q] 
+ riV2Q + r2(VV-Q)sT (5) 

to the lowest order in a series expansion in powers of VQ. 

Since Q is a traceless and symmetric second rank 3x3 
tensor it has five independent components. Accordingly, 
when the equation of motion of the alignment tensor are 
appropriately scaled, it is possible to express it in the 
following orthonormalized basis 

Q = y^QiTj, 

i 

To = v/372(zz)sT, 

Ti = v/i72(xi-yy), 

T2 = \/2(xy)sT, 
T3 = \/2(xz)sT, 

T4 - V2(yz)sT, (6) 

and study the equations of motion of each of the compo- 
nents Ofc, fc = 0, 1, . . . , 4||3^ projected out. 

It has been observed in the absence of spatial varia- 
tion that depending on the model parameters entering 
the equations, the order parameter equations can have 
different characteristic orbits [3 US- Possible in-plane 
states, where as the name suggests, the director is in 
the plane of flow determined by the direction of the flow 
and its gradient, and the order parameter components 
03, 04 = are "Tumbling" (T, in-plane tumbling of the 
alignment tensor), "Wagging" [W , in-plane wagging) and 
"Aligning" (A, in-plane flow alignment) states. Out of 
plane solutions, characterized by non-zero values of aa 
and 04, observed are "Kayaking-tumbling" {KT, a peri- 
odic orbit with the projection of the main director in the 
shear plane describing a tumbling motion), "Kayaking- 
wagging" [KW, a periodic orbit with the projection of 
the main director in the shear plane describing a wag- 
ging motion) and finally "Complex" (C) characterized 
by complicated motion of the alignment tensor. This in- 
cludes periodic orbits composed of sequences of KT and 
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KW motion and chaotic orbits characterized by a posi- 
tive largest Lyapunov exponent. 

A solution phase diagram based on the various in-plane 
and out-of-plane states ioi A — and ai = is given 
in[35i|. It is observed that ai ^ gives similar results [T^. 
As control parameters, we use Xk = — (2/V3)ao related 
to the tumbling coefficient in Leslie-Ericksen theorvfs^ 
Issf , and the shear rate 7 to study the phase behavior of 
this system. 

It is observed in experiments that the flow curve (shear 
stress vs strain rate) of a wormlike micellar system in 
shear flow has a rather large plateau region where banded 
flow is believed to occur, and a study of the dynamics 
of traceless symmetric order parameter Q (Eq^ for a 
sheared nematogenic system, that allows spatial varia- 
tion, is likely to capture this feature. As we shall see later 
the shear banding observed in such systems is dynamic in 
nature and is an important element in the spatiotemporal 
rheochaos we observe. 

Hereafter we express the equations in the orthonor- 
malised basis as in Eq|H| As in[33. Issf. we rescale time 
by the linearized relaxation time r/A* at the mean-field 
isotropic-nematic transition, and Q as well by its mag- 
nitude at that transition. We have set ai = in our 
analysis as it seems to have little effect on the dynami- 
cal behavior of the system^, Is^ . Issj in the parameter 
range studied. Further, we choose ^ = throughout, to 
make a correspondence to the ODE studies of |33. l35l| . 
This places the system well in the nematic phase at zero 
shear, in fact at the limit of metastability of the isotropic 
phase. Distances are non-dimensionalized by the diffu- 
sion length constructed out of Fi and r/A*. The ratio 
F2/F1 is therefore a free parameter which we have set to 
unity in our study. 

The resulting equations are then numerically inte- 
grated using a fourth order Runge-Kutta scheme with 
a fixed time step {At = 0.001). For all the results quoted 
here a symmetrized form of the finite difference scheme 
involving nearest neighbors is used to calculate the gra- 
dient terms. Thus 



(A.t)2 
fi+l — fi-1 



2Ax 



(7) 



We have checked that our results are not changed if 
smaller values of At are used. We have further checked 
that the results do not change if the grid spacing is 
changed (i.e. Ax is decreased) and more neighbors to the 
left and right of a particular site in question are used to 
calculate the derivative. This gives us confidence that the 
results quoted here do refiect the behavior of a continuum 
theory and are not artifacts of the numerical procedure 
used. We use boundary conditions with the director be- 
ing normal to the walls. With this, we discard the first 
6 X 10^ timesteps to avoid any possibly transient behav- 
ior. We monitor the time evolution of the system for the 
next 5 X 10^ time steps (i.e. t = 5000), recording con- 



figurations after every 10^ steps. We have carried out 
the study with system sizes ranging from L — 100 to 
L = 5000. 

Further, we calculate the contribution of the alignment 
tensor to the deviatoric stress [s^, 13 IHl c*^^ oc 
aoG — ai(Q • G)st where G, defined in (EqlSJ, is the 
nematic molecular field, and the total deviatoric stress is 
cr'^^ plus the bare viscous stress. Since the latter is a 
constant within the passive advection approximation, we 
can study the rheology by looking at cr'-'^ alone. We are 
aware of the importance of allowing the velocity profile 
to alter in response to the stresses produced by the order 
parameter field, and this is currently under study [s^. 

While generating the results for the time-series analy- 
sis for the LS, we run the simulation till t = 20000, for 
L — 5000 sized system, recording data at space points at 
intervals of I = 10. We monitor the space-time evolution 
of the shear stress (the xy component of the deviatoric 
stress cr'^^) (referred to as T,xy) and the first and sec- 
ond normal stress differences Yi^x — '^yy and Y.yy — Y,zz 
respectively. 



B. Results and Discussion 

1. Phase Behavior and Dynamic Shear Banding 

In view of prior work on observation of chaos in the lo- 
cal equations of motion of the alignment tensor (s patially 
homogeneous version of EqQ by Rienacker et al. |3J, |33 
we address the following question: Is the phase diagram 
(in the j — Xk plane) affected by allowing spatial varia- 
tion of the order parameter ? We answer this question in 
affirmative, and show that the 'C region of the phase di- 
agram of Rienacker et aLj^E^ corresponding to "com- 
plex" or chaotic orbits broadens upon incorporating the 
spatial degrees of freedom. In other words, there exist 
parameter ranges where the spatially homogeneous sys- 
tem is not chaotic, but chaos sets in once inhomogeneity 
is allowed. A result of particular interest is the obser- 
vation of "spatiotemporally intermittent" (STI) states in 
a certain range of parameters en route from the regular 
to the spatiotemporally chaotic regimes. Such behavior 
is by definition not accessible in the evolution equations 
of the spatially homogeneous alignment tensor studied in 
[33. l35| . It would be of interest to find chaotic regimes 
where only two of the five independent components of Q 
are nonzero. Since the number of degrees of freedom per 
space-point would then be two, such chaos would clearly 
be a consequence of spatial coupling. We have not lo- 
cated such a regime so far. 

Local phase portraits (orbits obtained when various 
pairs of order parameter components are plotted against 
each other) illustrate the chaotic or orderly nature of the 
on-site dynamics. Shown in the right panels of Fig|21are 
the local phase portraits (ai vs oq) for a particular point 
xq for various values of the tumbling parameter Afc, ob- 
tained by holding the shear-rate fixed at 7 = 3.5. We 



FIG. 2: Plots showing ao{x,to) vs ai{x,tQ) (left panel) 
and ao{xo,t) vs ai{xo,t) (right panel) for Periodic [(a),(e)], 
Chaotic [(b),(f)], (C^A) [(c),(g)] and aligned [(d),(h)] 
regimes. 

have checked that the character of the phase portrait re- 
mains intact upon going from one space point to another 
though in the chaotic regime there is no phase coherence 
between two such portraits. A closed curve correspond- 
ing to a limit cycle is seen at Xk = 1.27, while at Afe = 1.3 
(corresponding to the 'C region of the phase space) it is 
space filling. At — 1.35, as one approaches the region 
where the director aligns with the flow the points reduce 
to those on a line and eventually in the aligning regime 
(Afc = 1.365) where the director has already aligned with 
the flow it is represented by a point. This assures us 
that the local dynamics in the spatially extended case is 
similar to that of the ODEs of[3ll35t. 

We also construct the spatial analogues of these por- 
traits, i.e., we allow the system to evolve till a sufficiently 
long time (say to) and then record the spatial series. 
Again we get a limit cycle in the 'T' region of 0,123)) fol- 
lowed by a space-filling curve in the 'C region. As we go 
from the Chaotic towards the Aligning regime, the points 
arrange themselves on a line, and finally in the aligned 
regime one only obtains a point, corresponding to a spa- 
tially uniform state. This is shown in the left panel for 
FigEl One should note here that the 'T' regime here cor- 
responds to spatiotemporally periodic states whereas in 
other regions of parameter space one does observe states 
that are temporally periodic but spatially homogeneous 
(Fig. and for such states the local phase portrait 

corresponding to spatial variation at a fixed time would 
be a point and not be a closed curve. 

We now turn to the detailed spatiotemporal structure 
of the phase diagram of this system. We find many inter- 
esting phases including spatiotemporally chaotic states 
with a broad distribution of length scales (Figsn^e)|2J), 
spatiotemporally irregular states in which a few length 
scales are picked up by the system (Figs^f) and (g), 
0J, a flow aligned phase (Fig. H^h)) and also regular 
states (R) showing periodicity in both in time and space 
(Figl^l or that are periodic in time and homogeneous in 
space (Fig^a)). In addition to these states we find the 



FIG. 3: (Color online) Space-time behavior of the shear stress 
in the chaotic regime, 7=3.678 and Afc = 1.25. Slice taken from 
a system of size L = 5000. 




FIG. 4; (Color online) Space-time behavior (surface plots) 
of the shear stress in the chaotic to aligning regime, 7 = 4.05 
and Afc=1.25. Slice taken from a system of size L — 5000. 

presence of spatiotemporally intermittent states (STI) 
(Fig^b) and (c)). In regions of parameter space we have 
also observed spatially and temporally ordered domains 
co-existing with patches characteristic of STI.(Fig|Hl). 
The parameter values at which these are seen, corre- 
spond well with those obtained from the phase diagram 
of|3l|35|. 

Let us now focus on the parameter region labelled 'C 
or 'Complex' in|3^ l35l |. where we find spatiotemporal 
chaos. This regime is characterized by dynamic insta- 
bility of shear bands as seen in FigOl which shows sev- 
eral distinct events, such as the persistence, movement, 
and abrupt disappearance of shear bands. It is found 
that the typical length scale at which banding occurs is 
a fraction of the system size, though it follows a broad 
distribution. As one moves closer to the phase boundary 
separating the spatiotemporally chaotic state from sta- 
ble flow alignment, the bands become more persistent in 
time and larger in spatial extent as shown in Fig^ 

We have computed the distribution of band sizes or 
spatial "stress drops", and looked for the presence of 
dominant length-scales in the system in order to obtain 
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FIG. 5: (Color online) The space-time correlation function 
in a regime which is both spatially and temporally periodic. 
7=3.4,A=1.29 




FIG. 6: (Color online) Space-time evolution of the shear 
stress showing co-existence of different dynamical regimes at 
7=3.9,A=1.12. 

a better understanding of the disorderly structure of the 
shear bands as seen in Fig|31 and compared it with the 
behavior seen close to the phase boundary (Fig^. An- 
other important reason for such an analysis is to rule 
out any hidden periodicity that might be present in the 
space-time profiles of shear stress as shown in Fig 01 

The "stress drop" calculation is outlined below. At a 
given time (say U), we define a threshold Soxyj & little 
above the global mean {T^xy)y,t, and map the spatial con- 
figuration to a space-time array of ±1: T^xy = ngniT^xy ~ 
Soj;y). Figd shows the histogram of the spatial length 
of intervals corresponding to the -l-state, for the Chaotic 
and the Chaotic to Aligning (C A) regimes. We have 
considered configurations extending over L — 2500 spa- 
tial points, and the statistics is summed over configu- 
rations sampled at 5000 times (i.e. i — 1,5000). As ex- 
pected, the distribution of band lengths in the spatiotem- 
porally chaotic regime is fairly broad and roughly expo- 
nential in shape, whereas as one approaches the Aligning 
regime, the distribution is peaked about a few dominant 
length-scales. Also, note that as one passes from the 
Chaotic to the Chaotic-to-Aligning (C A) state, the 
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(c) , (d) 

FIG. 7: (Color online) Spatial distribution of "stress drops" 
(corresponding to residence intervals in which the shear stress 
is above a threshold Eo^y = 0.8) in the Chaotic (a) and C 
A (b,c,d) regimes. 7 = 4.0 and Afc = 1.22 (a), 1.24 (b),1.25 
(c), 1.27 (d). 

dominant length-scale associated with the shear bands 
increases. 



2. Route to the Spatiotemporally Chaotic State 

We now monitor the approach to the spatiotemporally 
chaotic state as a function of the tumbling parameter 
A/c, for a fixed value of 7 (= 3.8). We observe a se- 
quence of states. At low (1.0), the shear stress is 
periodic in time and homogeneous in space Fig^a). As 
we increase Afe, we come across states which are both 
spatially and temporally disordered, (Figs^lb) and (c)) 
consisting of propagating disturbances in a background 
of highly irregular local structures, which resemble ge- 
ometric patterns seen in probabilistic cellular automata 
[5^ . The borders of the ordered regions evolve like fronts 
towards each other until this region eventually disap- 
pears in the chaotic background. These states are typi- 
cal of what is known as "spatiotcmporal intermittency" 
(STI) 58, 59, 60, 62]. Indeed, it is suggested 60, 65] that 
the transition to fully developed spatiotcmporal chaos 
generally occurs via this admixture of complex irregular 
structures (high shear stress) intermittently present with 
more regular low shear regions. In contrast to low dimen- 
sional systems where intermittency is restricted to tern- 
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poral behavior, STI manifests itself as a sustained regime 
where coherent-regular and disordered-chaotic domains 
coexist and evolve in space and time. Earlier studies of 
the onset of spatiotemporal chaos ^5^, 'b^, "6^, |^ l6^ |6^ 
suggest a relation to directed percolation (DP). Evi- 
dence for DP-like behavior in spatiotemporal intermit- 
tency mostly comes from studies of coupled map lattices 
(CML) 62]. Such processes are modeled as a probabilis- 
tic cellular automaton with two states per site, inactive 
and active, corresponding respectively to the laminar and 
chaotic domains in the case of STI. Studies find that in 
the STI regime, a laminar (inactive) site becomes chaotic 
(active) at a particular time only if at least one of its 
neighbors was chaotic at an earlier time, there being no 
spontaneous creation of disordered-chaotic sites. Hence a 
disordered site can either relax spontaneously to its lam- 
inar state or contaminate its neighbors. This feature is 
analogous to directed percolation, and one consequence 
of this picture is the presence of an absorbing state: in 
STI studies of CMLs, once all the sites relax sponta- 
neously to the laminar state, the system gets trapped 
in this state forever, thus the laminar state in STI corre- 
sponds to an absorbing state in DP. This analogy predicts 
that STI should show critical behavior similar to that as- 
sociated with DP - power law growth of chaotic domains, 
and characteristic static and spreading exponents. There 
is however still no uniformity of opinion on whether spa- 
tiotemporal intermittency belongs to the same universal- 
ity class as DP as characterised by the critical exponents 
of the DP class. Some studies of coupled map lattices 
[58l |. PDEs |6^, and experiments j6^ suggest that, 
though the critical behavior in STI is visually similar to 
DP, the exponents measured in STI are not universal. 
Other investigations that have evaluated the exponents 
at the onset of spatiotemporal intermittency in coupled 
circle map lattices and in experimental systems |6^. 
claim that this transition indeed falls in the universal- 
ity class of directed percolation. Our aim in this context 
is to macroscopically characterize the disordered struc- 
tures in the spatioteniporally intermittent state in terms 
of the distribution of the widths of the laminar domains 
and study the qualitative connection between STI and 
DP. Accordingly, following the analysis in Chate et al. 
[ssl l65l | we study the decay of the distribution of lami- 
nar domains at the onset of intermittency and far away 
from this onset. At each time-step, the spatial series of 
the shear stress values are scanned, and the widths of 
the laminar regions (regions for which the spatial gradi- 
ent is less than a sufficiently small value) are measured 
and inserted into a histogram. This process is then cu- 
mulated over time, giving the distribution of laminar do- 
mains. Previous studies have found that at the onset 
of spatiotemporal intermittency, this distribution has a 
power law decay (with a power ranging from 1.5 to 2.0 
in CML studies [Sg] and experiments 1691 and 3.15 in a 
variant of the Swift- Hohenberg eqn. |65j|'). while away 
from the onset, it has an exponential decay. We find a 
similar behavior. In our work, at a representative point 




FIG. 8: Decay of the distribution of laminar domains at the 
onset of the spatiotemporal intermittent regime (represented 
by circles) on Log-Log scale. The solid line is a power law fit 
to the data with an exponent tpsTi = 1.86 ± 0.05. 




150 



FIG. 9: Decay of the distribution of laminar domains away 
from the above-mentioned onset and well within the STI 
regime (represented by circles) on a semiLog scale. 



(7 — 3.9, Afc = 1.116) at the onset of intermittency, the 
distribution of laminar domains has a power law decay 
(Fig. inj, with an exponent of 1.86 and a standard de- 
viation from the data of 0.05. Away from this onset 
(7 — 4.0,Afc = 1-13), the decay is close to exponential 
as evident from Figure El 

Another way of characterizing the time evolution of the 
coherent and disorderly regions in the spatiotemporally 
intermittent states of Figsn(b) and (c) is to calculate the 
dynamical structure factor for the shear stress. FigllOl 
shows the dynamic structure factor characteristic of this 
regime. When the system is in the fiow aligned regime 
S{k, Lo) has a peak at fc = and w = 0. On the other hand 
in the spatiotemporally intermittent state, the dominant 
weight in S(k,uj) is on lines of w oc ±fc, implying distur- 
bances with a characteristic speed of propagation. The 
front velocity of the cellular automata like patterns seen 
in the space-time plots in this regime can be calculated 
from the slope of these lines. In FigsQ (d) and (e), the 
system is in the chaotic regime. As we pass on from the 
chaotic towards the aligning regime, more regular struc- 
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FIG. 10: (Color online) Pseudocolor plot ol the dynamic 
structure factor S{ky,u!) in the "STI" (spatiotemporally in- 
termittent) regime. Color at any point ^ logarithm (to base 
10) of the value of S{ky,uj) at that point. 
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FIG. 11: (Color online) Dynamic structure factor in the spa- 
tiotemporally periodic regime. The peridicity in time and 
space is borne out by the straight lines parallel to the fre- 
quency and wave-vector axes respectively. 



tures are seen to evolve (Figsn(f), (g)); the shear bands 
grow in spatial extent and are more long-hved. Fig^ (h) 
shows a snapshot of the shear stress in the flow aligned 
regime. Figures ITTl and IT^ show the dynamical structure 
factor in the spatiotemporally periodic regime and close 
to the chaotic-aligning phase boundary respectively. 

Finally, we present the phase diagram coming out of 
the previous analysis as in the FigEl The phase dia- 
gram was calculated by monitoring the shear stress pro- 
file as well as various components of the nematic order 
parameter (as they are strikingly different for the differ- 
ent phases, which is also borne out by the difference in 
the corresponding dynamic structure factors) as a func- 
tion of the shear rate 7 and the tumbling parameter Afc. 
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FIG. 12: (Color online) Dynamic structure factor in the C 
A regime. The system has almost relaxed to a steady state 
in time, and spatially there are large domains that are flow 
aligned. 
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FIG. 13: (Color online) Phase diagram of the system in the 
\k vs 7 plane, showing regular, i.e., periodic in time and either 
periodic or homogeneous in space (R), chaotic (C), aligning 
(A) and spatiotemporally intermittent (STI) regimes. Note 
the reentrant chaotic behavior as a function of 7 in a narrow 
region of and as a function of for 7 < 3.6. 



3. Lyapunov structure of the Chaotie State 

Next we try to characterize the chaotic states in our 
study. In studying dynamics of spatio-temporal systems 
[68l | , one needs to establish whether the system is truly in 
a spatiotemporally chaotic regime or can be described by 
a model with only a few (dominant) independent modes. 
So from the multivariate time-series generated by such 
systems, one tries to compute quantities analogous to the 
invariant measures used to characterize low dimensional 
chaos. However, true spatiotemporal chaos corresponds 
to spatially high dimensional attractors, with dimension 
growing with the system's spatial extent, and the esti- 
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mation of invariants such as the correlation dimension 
can be quite problematic. Indeed we find that the chaos 
that we observe is quite high dimensional (embedding 
dimension '6^, TO > 10). A reliable estimate of the 
correlation dimensions can be made only from a data 
train so long as to require prohibitively large computa- 
tional times to generate. 

An alternative approach is to study the Lyapunov spec- 
trum (LS). For a discrete N dimensional dynamical sys- 
tem, there exist N Lyapunov exponents corresponding 
to the rates of expansion and/or contraction of nearby 
orbits in the tangent space in each dimension. The LS 
is then the collection of all the A'' Lyapunov exponents 
Xi,i = 1 : N, arranged in decreasing order. The LS is 
very useful in the characterization of a chaotic attrac- 
tor. Useful quantities that can be calculated from the LS 
are the number of positive Lyapunov exponents Nx^ and 
sum of the positive Lyapunov exponents X]a+ ■ ^^'^^ 
the sum of the positive Lyapunov exponents provides an 
upper bound for the so called Kolmogorov-Sinai entropy 
h which quantifies the mean rate of growth of uncertainty 
in a system subjected to small perturbations. In many 
cases, h is well approximated by the sum J2\ 71]. Both 
these quantities have been found to scale extensively with 
system size in spatiotemporally chaotic systems. For dy- 
namical systems with only a few effective degrees of free- 
dom, it is straightforward to compute the LS. However 
for extended systems with a large number of degrees of 
freedom, even a few hundred, it runs into severe diffi- 
culties because of the inordinately large computing time 
and memory space required. In such situations it is im- 
portant to make use of techniques that derive informa- 
tion about the whole system by analyzing comparatively 
sma ll sy stems with exactly the same dynamical behavior 
(zIIzIItII. It has been widely observed that the LS for 
spatiotemporal systems is an extensive measure [t^ and 
is associated with a rescaling property [z3ilzlil3 i-e., the 
LS of a subsystem, when suitably rescaled can give rise 
to the LS of the whole system jTg • The volume rescaling 
property for the LS in spatiotemporally chaotic systems 
also implies that extensive (size dependent) quantities 
such as X]a+ ^^'^ scale with not only the system 
size but also the subsystem size. Hence, instead of trying 
to study the spectrum and related quantities in a system 
of large size N, one could confine the analysis to rela- 
tively small, more manageable subsystems of size i.e., 
at space points j in an interval io < j < io+Ns — 1 (where 
io is an arbitrary reference point), and study the scaling 
of related quantities with subsystem size A^sl6£|. Thus, 
instead of trying to implement the correlation-dimension 
method for our spatially extended problem, we study the 
LS [1^123 ■ Further, instead of studying systems of ever- 
increasing size, we look at subsystems of size Ng in a 
given large system of size A^. 

For spatiotemporal chaos we expect to find that the 
number of positive Lyapunov exponents grows systemat- 
ically with Nfi- This is seen in Fig^] For both figures 
in Fig ll4l , we carry out the procedure with two different 
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FIG. 14: Sum of positive Lyapunov exponents (left panel). 
Number of positive Lyapunov exponents (middle panel) and 
the largest Lyapunov exponent (right panel) as functions of 
subsystem size Ns, for 7 — 3.678, Xk ~ 1.25. Embedding 
dimension for the time series of each space point is 10.(io=101, 
see text). 



reference points io and find essentially the same curves. 
Furthermore, it has been reported in many studies of 
spatiotemporally chaotic systems [t^, that when 

calculating the subsystem LS for increasing subsystem 
size A^s, one finds that the Lyapunov exponents of two 
consecutive sizes are interleaved, i.e. the zth Lyapunov 
exponent for the sub-system of size A^, lies between the 
ith and (z -I- l)th Lyapunov exponent of the subsystem of 
size A^s -|- 1. A direct consequence of this property is that 
with increasing subsystem size A^^, the largest Lyapunov 
exponent will also increase, asymptotically approaching 
its value corresponding to the case when the subsystem 
size is of the order of the system size. This trend is clearly 
seen in Fig^^ (b). 



C. Conclusions 

In summary, we have proposed a mechanism by which 
one might explain the chaotic and irregular rheological 
response of soft materials in shear flow, wormlike mi- 
celles in particular. The main idea brought out in this 
paper is that the coupling of orientational degrees of free- 
dom in a complex fluid with hydrodynamic flow can lead 
to spatiotemporal chaos for low Reynolds number flows. 
In particular, we have demonstrated that the nonlinear 
relaxation of the order parameter in nematogenic fluids, 
together with the coupling of nematic order parameter 
to flow, are key ingredients for rheological chaos. The 
broad idea that nonlinearities in the stress and spatial 
inhomogeneity are essential is a feature that our work 
shares with ^ and ;37t i38i] . 

We should note here that there could be more than one 
mechanism at work in producing rheochaos. The mech- 
anism observed in our study is that the system exhibits 
chaos in its local temporal dynamics, and then these lo- 
calized regions mutually interact with one another to gen- 
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erate spatial disorder. Fielding et al. 9] however have 
found spatiotemporally chaotic rheological behavior in a 
model system whose local dynamics does not show chaos, 
but incorporating the spatial degrees of freedom makes it 
chaotic. Note also that in the parameter range we have 
studied so far the equilibrium phase of the system, in the 
absence of shear, is neniatic. We made this choice to facil- 
itate comparison with the work of ^34. ,35] ; a better choice 
from the point of view of the experiments on wormlike 
micelles would be to work in the isotropic phase, with a 
substantial susceptibility to nematic ordering. We do not 
know if shear produces chaos in that situation, although 
it seems likely. It is also worth investigating whether ne- 
matics with stable flow- alignment at low shear-rates can 
go chaotic at higher rates of flow. That the structures in 
the transitional region between order and spatiotempo- 
ral chaos are similar to those in directed percolation, as 
in some other systems undergoing the transition to spa- 
tiotemporal chaos, is interesting and suggests a direction 
for possible experimental tests. 

We now comment on experiments which can test some 
of the ideas proposed in this paper. The dynamics of 
the alignment tensor can be studied in rheo-optical ex- 
periments on dichroism|78j , flow birefringence and rheo- 
small angle light scattering [f^. Flow birefringence ex- 
periments carried out in the last decade have shed light 



on shear banding and orientational properties of micel- 
lar solutions 0- Small angle neutron scattering experi- 
ments, using a two-dimensional detector, have also been 
used to analyse the orientational degrees of a micellar 
fluid in shear flow: the presence and proportions of the 
isotropic and nematic phases under shear, as well as the 
order parameter of the shear induced nematic phase in 
such systems have been studied In order to investi- 
gate rheochaotic behavior in space and time in systems 
such as under consideration, one could use these rheo- 
optical techniques and try to look for the irregularities 
in the spatial distribution of band sizes and their tempo- 
ral persistence, in a regime in the nematic phase where 
the micelles are not flow aligned. Further, very recently, 
spatio-temporal dynamics of wormlike micelles in shear 
flow has been studied using high-frequency ultrasonic 
velocimetry.23:], and various dynamical regimes includ- 
ing slow nucleation and growth of a high-shear band and 
fast oscillations of the band position have been observed, 
though the complex fast behavior reported is not chaotic. 

We thank G. Ananthakrishna and R. Pandit for very 
useful discussions, and SERC, IISc for computational fa- 
cilities. MD acknowledges support from CSIR, India, and 
CD and SR from DST, India through the Centre for Con- 
densed Matter Theory. 
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